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Abstract. We introduce a chipping model in one dimensional periodic lattice with 
continuous mass. The model evolves following a parallel or random sequential dynamics 
where a fixed fraction of the mass of any site is distributed randomly among the 
departure site and its neighbors, and the remaining mass sticks to the site. We have 
calculated the steady state mass distribution of the model perturbatively for both 
symmetric and asymmetric mass-transfer. In most cases, the product measure turns 
out to be a good approximation. 
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1. Introduction : 

Most systems in nature are in non-equilibrium states [1], in a way that the accompanying 
fluxes of mass, energy, or spin etc. are irreversible. Unhke their equihbrium counterparts 
where the stationary state is characterized by the Gibbs measure, these systems usually 
reach different and novel stationary states depending on the dynamics of the microscopic 
constituents. Several non-equilibrium lattice models have been proposed recently [21 13] 
to investigate the unusual steady state distributions, spatio-temporal correlations and 
possibility of macroscopic collective phenomena. 

One of the simplest non-equilibrium model is the mass transfer model where each 
site of a lattice is associated with discrete masses (particles) following a dynamics that 
involve aggregation, fragmentation, adsorption or desorption [U El El [TJ [8] . Interestingly, 
these systems undergo a condensation phase transition as the density of the system is 
increased. Study of these model systems have generated considerable interest among 
physicists, as a wide variety of natural systems exhibit basic microscopic mechanism 
similar to that of the simple mass transport models. These includes colloidal suspensions 
[9], polymer gels [lOl [H], river networks [12], traffic models [l3], wealth distributions 
[E], abelian sandpile models of self organized criticality [15], zero range process [I6] etc. 

Recently, a continuous version of the mass transfer model has been proposed by 
Rajesh et al. [IT] and some of its variations has been studied by others [iHlllQj. Unhke 
the discrete mass model, where a fixed unit of mass is transferred to the neighboring sites, 
here a random fraction of the mass is allowed to be redistributed among the sites and 
their neighbors [17J. Many of these models are known to evolve into a non-equilibrium 
steady state that is factorized. 

In this article we introduce stickiness, quantified by a parameter A, to the continuous 
mass transfer models. At each site, (1 — A) fraction of the mass is chipped off (thus A 
fraction of the mass sticks to the site), which is then redistributed either asymmetrically, 
i.e. among the departure site and its right neighbor, or symmetrically, i.e. among the 
two neighbors. Accordingly, the model is referred to as asymmetric or symmetric sticky 
chipping model respectively. We use a novel perturbation approach to calculate the 
steady state mass distribution of these models, for both parallel and random-sequential 
dynamics. Although the perturbation approach completely relies on the assumption that 
the steady state has a product measure, the resulting mass distributions are strikingly 
close to that obtained from Monte-Carlo simulations. 

The article is organized as follows. The model and the perturbation method is 
outlined in section [2l In section |3] we study the asymmetric version of the model and 
obtain the steady state distribution. The symmetric version of the model, where the 
chipped off mass is distributed among both the neighbors, is discussed in section HI 
Finally we discuss the results in section El and conclude. 
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2. The Model : 

The model is defined on a one dimensional periodic lattice with sites labeled by 
i = 1,2 . . . L. A continuous variable Xj, called mass, is associated at each site i which 
evolves according to the following dynamics. In each time step, (1 — A) fraction of 
the mass is chipped off from a site i (thus, A fraction of the mass sticks to the site) 
and then it is distributed among the departure site i and its neighbors i ± 1. In this 
article we study two different dynamics, (a) where the chipped off mass (1 — X)xi is 
randomly distributed among the two sites, i and its right neighbor i + 1 and (b) where it 
is distributed randomly among the neighbors {i + 1) and {i — 1). To distribute the mass 
randomly among two sites we draw a random number r from a uniform distribution and 
then allow one of the sites to take r fraction of the chipped mass; the rest, i.e. (1 — r) 
fraction of the chipped mass, is given to the other site. Hence forth we refer to model (a) 
as the asymmetric sticky chipping model (ASCM) as the mass is always transported to 
the right, whereas model (b) will be referred to as the symmetric sticky chipping model 
(SSCM). Both these models are studied with parallel and random sequential dynamics. 

These models differ from the usual chipping models as for non-zero A, there is a 
upper bound on what amount of the mass can be chipped off from a site, as each lattice 
site is constrained to retain a fixed fraction of the available mass. We must mention 
that for A = 0, in model model (a), r fraction of the whole mass Xi is given to the right 
neighbor. This model is same as asymmetric chipping model studied earlier [17], for 
both parallel and random sequential dynamics. 

For ASCM and SSCM, we could calculate all the moments for both parallel and 
continuous time dynamics assuming that the steady state has a product measure, which 
is known to be an excellent approximation (though not exact) [19j in many models. The 
closede form expression for the steady state mass distribution P(x) could be obtained, 
for A = in all cases, and for A 7^ in a few specific cases. 

In our effort to construct P{x) from the moments, we use a novel perturbation 
approach. The general principle of this approach is described in the following subsection. 

2.1. Perturbation approach I : 

Here, we assume that all the moments of the steady mass distribution are known 
and then try to construct P{x) perturbatively. Although P{x) depends on A, we omitted 
the argument for notational convenience. Let us express P{x) as a power series in A, 

about A = where we can obtain the exact expression of P{x), 

00 

P{x) = Poix) + \P^{x) + \^P2{x) + . . . = E >^'Pk{x) (1) 

A:=0 

where functions Pk{x) do not depend on A. In the chipping models discussed here, the 
total mass of the system and hence the density (x) = j^Y^iXi is conserved. Without any 
loss of generality one can fix the density to be unity. This imposes a condition on P{x) 

/■oo /"OO 

/ dx xPix) = (x) = 1= dx Pix) , (2) 

JO JO 
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where the last equahty stands for the normahzation condition. Now, for A 
P{x) = Po{x). Therefore, Po{x) satisfies two conditions, 

/ dx Pq{x) = 1 and / dx xPq{x) = 1. 

JO JO 

Thus for any other /c > 0, using (jj]) and ([3]) we have 



dx Phix) 







dx xPk{x). 



(3) 



(4) 



Since above constraints can not be satisfied by a real positive function (as x > 0) one 
can not interpret Pk{ probability distribution function (PDF). However, one can 

then obtain P(x) using Eq. ([1]), which is assured to be a strictly, positive function, even 
though Pk{x) are not. Pk{x) can be obtained directly from knowing the moments of 
P{x). To calculate Pk{x) let us expand the moments (x") as a power series in A, 



n QQ OO 

(x") = / dx x"P(x) = V Cjrh''. 

•^0 k=0 



(5) 



Here ^ are constant co-efficients (do not depend on A), which can be determined from 
Eqs. O and (ED, 



in) 



dx x^Pkix). 



(6) 



As discussed above, Pk{x) for any /c > is not a positive function, and thus should 
not be regarded as its n*^ moment. Pk{x) can be calculated from Cj^^ as follows. 



P{x) = 



Ln=0 



and from Eqs. ([7]), ([5]) and ([T]) it is evident that 



Ln.=0 



(7) 



(8) 



where =Sf ^ denotes the inverse Laplace transform. 



2.2. Perturbation approach II : 

Here we assume that the mass distribution P(x) satisfies a differential or an integral 
equation; the explicit form of the moments of the distribution P(x) are not known. In 
this case, one can first use the Laplace transform =Sf on the integral equation, which 
usually translates into a differential or a transcendental equation in (5(s) which is the 
Laplace transform of P(x), 

/■oo 

Q(s) = 5e \P{x)\ = / dx e-'^P(x). (9) 

JO 

We proceed further by expanding Q{s) as power series in A about A = 0, 

oo 

Q{s) = J2\'Qk{s), (10) 

fc=0 



Conserved mass models with stickiness and chipping 



5 



and equate the coefficients of different powers of A order by order. Finally one can find 
the distribution 



Pix) = \Q{s)\ = ^ 



-1 



.k=0 

These two approaches are equivalent as it is evident from Eq. ([8]) 



Q,is)=^[P,{x)] = j:'—^C, 



(n) 



ra=0 



(12) 



For any particular problem we will use the approach which is convenient. For AS CM, we 
show the equivalence explicitly, whereas for the SSCM we describe one of the approaches 
in details. 



3. Asymmetric Sticky Chipping Model (ASCM) : 

In this section we study the asymmetric version of the model, namely ASCM. Here, 
in each time step, A fraction of the available mass Xi at lattice site i remains sticked to 
the site. From the rest (1 — X)xi, fraction is transported to the right neighbor i + 1, 
and (1 — Ti) fraction is retained at the departure site i. In other words (1 — A)rjXj is 
transported form site i to i + 1. Here is a random number uniformly distributed in 
the interval (0, 1). 

First let us consider the parallel dynamics, where all the sites are updated 
synchronously. The dynamics can be written as 

Xi{t + 1) = \xi{t) + (1 - A)(l - r,)xi{t) + (1 - A)ri_ia;,_i(t) (13) 

for all 2 = 1,2 ... L, where the ffist term represents the mass that sticks to the site i, 
the second term corresponds to the mass that retains at the site i after (1 — X)riXi{t) is 
transported to i + 1. The third term results from the mass that site i receives from site 
i - 1. 

A special case of the model, with A = has been studied earlier [17] for both parallel 
and random sequential dynamics. It was shown that only for the parallel dynamics the 
steady state has a product measure Ps{{xi}) = Yli=i P{xi). For random sequential 
dynamics, a factorized steady state was only an excellent approximation. 

In this study, we assume that Ps{{xi}) is factorized, and try to obtain P{x) using 
one of the perturbative approaches discussed in section [2l The results are then compared 
with the stationary state mass distribution obtained from the Monte-Carlo simulations 
of the model. 



3.1. Approach I : 

In order to find the steady state mass distribution, ffist we calculate its moments. 
Since, in the steady state, distribution of x(t + 1) is same as the distribution of x{t), 
one may obtain (x") using Eq. ( JT3l) as follows, 

(^") = t Q ((A + (1 - A)(i - n)y){x'Ki - (rtt) (^r-t) . 
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This can be simplified further 



(n + iy. 



n 



(n+ 1)A + A"+^ - ( 
^(1- A'^+^)(l- A)"-'= 



n-l 



1 - A)"+i 



-1 



,=1 {k + l)\in-k + l)\- ^^^^ 

Since, the n*^ moment (x") here depends on all the lower moments up to n — 1, we 
calculate (x") recursively starting from (x) = 1. The first few of them are, 

3(A + 1) 



4A + 2 

3(A^ + 3A + 2) 
2(2A + 1)2 
15(A + 1)2(A2 



A + 3) 



2(2A+ 1)2(2A3- A2 + 6A + 3) 



(15) 



The moments (x*^) as a function of A becomes messy with increase of n, and obtaining 
a general expression becomes practically impossible. Thus it would be useful to obtain 
P(x) perturbatively which gives all the moments correctly up to some order in A. 

To proceed with the perturbative approach, first let us express (x"") as a power 
series in A as done in Eq. ([5]), 



A:=0 k=0 

where we have used, for convenience 



k ^ 1 



c, 



in) 



(16) 



(17) 



From the Taylor's series expansion of Eq. f ll4p in A about A = 0, we equate the 
coefficients of different powers of A and get a set of recursive equations for A^^\ It 
is evident that A^""* depends on A^^ ■* with < k' < k and < n' < n. These equations 
can be solved recursively as they are supplemented by the boundary conditions 



A 



(0) 



(0) 



A 



(1) 



2 ° 



where C, 



(0) 



'0 

Cq^^ corresponds to the normalization of the steady distribution 



function and the conservation of mass respectively (see Eqs. ([3]) and ([6])). We can then 
obtain 



Pk(x) = if- 



E(-^)"(r^ + l)A 



in) 



.n=0 



(19) 



using Eq. (|8]) and finally calculate P(x) using Eq. (IT| 
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Qth Qj^dg^ perturbation in X : 

In the perturbative approach the O*'* order corresponds to the system with no 
stickiness. In this case, Eqs. ( IT4|) and ( !T6|) results in 

n-1 

{n - 1)4"^ = J2 V n > 2, (20) 

k=l 

which needed to be solved using the boundary condition A'q'^ = | (which is equivalent 
to fixing the density of the system to unity). This results in 

4") = 1 V n > 0. (21) 

Further using Eq. ( |T9l) we get the steady state distribution function as 

Po{x) = Axe-'^\ (22) 

The mass distribution function Eq. ( !22l) is identical to the steady state distribution 
obtained earlier for the model with no stickiness {i.e. for A = 0) where the steady state 
is known to have product measure [T7] . 

In Fig. [T]we show P{x) obtained from Monte-Carlo simulation of a system of size 
L = 100 for different A. Fig. W(Jo) compares P{x) for A = with Eq. (!22l) (solid line). 




Figure 1. (a) Steady state mass distribution P{x) of ASCM with parallel dynamics 
obtained from Monte-Carlo simulations of system size L — 100 for A = 0, 0.1, and 0.2. 
(6) P{x) for A = (dashed line) is compared with Eq. (|22p . 

In presence of stickiness, we will show that P{x) can be obtained with reasonable 
accuracy when perturbative correction terms are added to Po{x) order by order. 



1"** order perturbation in A ; 



Next we proceed to calculate the distribution function up to the first order in 
A. Using Eq. (fT6|) and comparing the coefficients of A in Eq. (fT^ we obtain 

n-1 

{n - l)Ai") = ^[24')4"-') -{n- fc)4')4""')] 
k=i 
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which are needed to be solved using boundary conditions = 0, A\ = (obtained from 
Eqs. (jl]), and ( 1T7|) ). Solving the above equation, we obtain 

4") = 2-"n(l - 7 - + 1)) 

where 7 is the Euler constant, 

7 = 0.57721... 

and = ¥g is the Digamma function [201 [22] • Using Eq. (IT^ we get 



PUx) =££-^ 



.n=0 



-CI" 



^-1 



-2a; 



-X \^g{x) — 1} + xg{x) — 



,n=0 
1 



(n) 



4 



(23) 
(24) 



(25) 
(26) 



where g{x) = ln(2x) + 7 

and 7 is given in Eq. (I24p . 

To check how good is the 1'^* order correction, we compare P{x) = Po{x) + XPi{x) 
with the steady state distribution obtained from the Monte-Carlo simulations of ASCM 
on a lattice of size L = 100, for two different A = 0.1, 0.2 (see Fig. [21(a) and (b)). 
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Figure 2. (a) and (b): The steady state mass distribution obtained from Monte-Carlo 
simulations (dashed line) of ASCM on a lattice of size L = 100 with parallel dynamics 
for A = 0.1,0.2 are compared with the 1** order perturbation results (solid lines), (c) 
and (d) : Comparison of the same by the taking 2"^^ order correction. 



Clearly for A = 0.1, the perturbative result is in good agreement with P{x) obtained from 
the Monte-Carlo simulations, despite the fact that we have assumed product measure 
for the steady state. For A = 0.2, Po{x) + XPi{x) shows a small departure, which can 
possibly be substantiated by adding the 2"*^ order correction. 
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To calculate the 2 order correction term in A we proceed in a similar way; first 
we obtain the recursion relation for A.*-"^ 



2 ) 



n-l 



{n - 1)4"^ = 



k=l 



02 11 V /Ol 



.(1) .(n~l) I n{n + l) r^) 



(27) 

Since A^^^ and A^^^ are known function of n (see Eqs. (12T|) and (l23ll ). the above equation 
can be solved by using 



^2(.) = E(-^r4"\ 

n=0 



(2J 



which is the generating function for A'^\ Substituting for Ag"'' and A^""^ in Eq. ([2] 
multiplying both sides of the equation by (— s)" and taking the sum over n we obtain a 
differential equation for V2(s); 

s{s + 2YV^{s) + (s - 2)(s + 2)V2(s) 



+ 2s + 4 In^ (^) + 2(s - 4) In (^) + 



2 2 ' ^ (29) 

This equation can be solved using the boundary condition V^'(O) = (obtained from the 
condition (x)=l), resulting 



V2{S) 



{2 + s) 



s^ + As- 2s In (l + -) + (s - 2) In^ (l + 



(30) 



Now ^2"'* is the coefficient of (— s)" in the power series expansion of V2{s) 



At' 



(n+i)2"+i L " 2 ~ + ""'^^ + + 2''' + 2^'{-2 + 7 + (-1 + 7)^}^^'^ 



where Hi^^ are Harmonic numbers 1201. Thus, 



2e 



-1 



E(-^)"(n + i)4' 



(n) 



.n=0 



-2x 



4x^(/(x) - 2g{x) -\ + ^)- 2a;2(5/(x) - 4^(x) -5^ + 8 



+2x(2^2(^) ^ 3^^^) _ _ {^(3;) + 1} 



+ 25{x) 



(31) 



where g{x) = ln(2x) +7. Note that in the above expression the term 2S{x) is essential to 
assure that /q°° dx P2{x) = 0. In Fig. [2](c) and (d) we compare ^0(2^) + APi(x) + A^P2(a^) 
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with P{x) obtained from the Monte-Carlo simulations of the dynamics fll3p on a ring 
of size L = 100 for two different A = 0.1, 0.2. It is evident that for large A the 2"'^ order 
correction provides a substantial improvement over 1'^* order correction. 



3.2. Approach II : 



The steady state distribution can also be obtained from the second approach 
discussed in the subsection 12.21 To obtain an equation that P{x) must satisfy in steady 
state, we proceed by writing a master equation, for the dynamics ( |T3l) . 
f^P{x'i} = -P{x'^) + /o°° dxi^i /o°° dXi /o dri Sq dri_if{ri)f{ri_i)P{xi, 

X 6[x[ - {A + (1 - A)(l - r,)}x, - (1 - A)r,_ix,_i), 
where P{xi,Xi-i) is the two point joint distribution function. In the steady state. 



P(x') 



dx 



i-l 







dXi / dri / dri_if{ri)f{ri_i)P{xi,Xi^i) 



X Si^x'i - {A + (1 - A)(l - ri)}xi - (1 - A)ri_iXi_i) . 

For A = 0, it has been already indicated in [1^ that both the parallel and random 
sequential dynamics of the model can be treated by taking the probability distribution 
function of as 

f{r,)=p+{l-p)S{n) (32) 

where < p < 1 is a fixed number. Effectively the random variable is chosen to be 
with probability (1 — p) or otherwise it is chosen uniformly from the interval (0, 1). 
For p = 1 the dynamics corresponds to the parallel update and in the p — ?• limit the 
dynamics is equivalent to the random sequential update. 

Let us assume that the steady state mass distribution of the model has a product 
measure, i.e. there is no correlation among the sites. Thus 

roo roo rl rl 

P{x'i) = / dxi^i / dxi / dri / dri_if{r.^f{ri_i)P{x.^P{xi^i) 
Jo Jo Jo Jo 

X 6[x', - {A + (1 - A)(l - ri)}x^ - (1 - A)r,_ix,_i). (33) 
Now, taking the Laplace transform on both sides of the above equation we obtain 
Q{s) [~p{l- p)v{sil - A)) - (1 - pY + 1 

= pV(s(l - A)) +p(l -p)] y rfr Q(sr + As(l -r)) (34) 

where Q{s) is the Laplace transform of P{x), and V{s) is defined by 

d 



V{s) = I dr Q{sr), or Q{s) = —{sV{s)} 
Jo ds 



(35) 

as 

In fact Q{s) can be obtained from Eq. flMl) . and one can rewrite it only in terms of the 
function V as 



d 
ds 



{sV{s)}[ - p{l - p)v(s{l - A)) -{l-pf + l 
pV(s(l - A)) + p{l - p)] f dr £ ^^'^^ " ""^^^ 



m=0 



ml 



ds"^+^ 



{sV{s)} 
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where we have expanded Q(^sr + Xs{l — r)^ in Taylor series about A = 0, and used Eq. 
fl35|) . Since this differential equation incorporates infinite number of terms, it is not easy 
to find a general solution. We again use the perturbative approach as follows. If the 
argument of the function V depends on A, we first get a Taylor series expansion of V in 
A about A = 0, 

1^(^(1 - A)) = E ^l^m* (3T) 

Each of the V{s) and their derivatives are then expanded using 

oo 

y(s) = ^AVfc(s). (38) 

k=0 

Now using Eqs. (1571) and (155]) in Eq. (15^ and collecting the coefficients of different 
powers of A, order by order, one can obtain a set of differential equations in terms of 
Vk{s). 

Finally we proceed to calculate Vk{s) for specific dynamics, i.e. parallel {p = 1) 
and random sequential {p — i- 0) by using the boundary conditions, 

k = : \/o(0) = l , K,'(0) = -^ 

k^O : Vk{0) = , V^{0) = (39) 

which comes from the generic conditions (jS]) and (jlj), correspond to the normalization 
of the probability distribution function and that the mass density is unity. 
Some comments are in order. Note that from Eq. ( l35i) we have 

Vk{s) = f' dr Q,{sr) or Q,{s) = ^{sVk{s)}, (40) 
JO ds 

where Q(s) = Ek=o ^"Qkis). Again, since Qk{s) = E^=o(-s)"(^ + 1)4"^ 

oo 

Vk{s) = j:i-^r^k^- (41) 

n=0 

Thus Vfc(s) is simply a generalization of V2(s) defined in Eq. (1281) . In other words Vfc(s) 
is the generating function of A^j^\ 

Now, let us consider the specific case p = 1 corresponding to parallel dynamics and 
proceed from Eq. fl36l) to obtain the coefficients of different powers of A, order by order. 

Qth Qj-^Qj. perturbation in A ; 

The terms independent of A in Eq. fl36l) gives 
- K)(s)' + Vo{s) = 
which is solved using the boundary condition 



V.(s) = (42) 
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Po{x) = 



4X6 



-2x 



Po{x) is identical to the exact mass distribution function obtained earher [TT] for A = 0. 

1^* order perturbation in A ; 

Equating coefficient of A in Eq. (1361) . we get 

s{Vo{s) - l}V^{s) + sVlis) + {1 - 2Vo{s)}V,{s) = 
where Vo(s) is given in Eq. fH2]) . This can be solved using boundary condition (15^ . 

2.1n(^) 



which results in 



PUx) = 



s + 2)2 
d 



(43) 



ds 



-2x 



—x'^{g{x) — 1} + xg{x) 



(44) 



where g{x) = ln(2a;) + 7 and 7 is Euler constant. Thus we obtain the same 1** order 
correction as in Eq. fl25|) obtained using approach 13.11 



2 order perturbation in A .■ 

The coefficients of A^ in Eq. gives 

-\v^{s) [s {sV^\s) - 2Vi{s) + 1} + AV^is)] + sV^{s)V^{s) - sV^is) 

~sVl{s) - Vi{sf + sV^{s) + V2{s) = 

where Vo(s) and Vi(s) are given by Eq. P2l) and (143|) respectively. One can solve this 
differential equation subjected to the boundary condition (j39l) . 



(2 + s)3 

Thereafter P2{x) is obtained as 

d 



+ 4s - 2s In 1 + - + (s - 2) InM 1 + - 



P.ix) = 



ds 



2e 



-2x 



AxHg'^ix) - 2g{x) - ^ + 3) - 2xHbg\x) - 4^(x) - 5^ + 



TT 



6 



6 



+2x{2g\x) + ?,g{x) - 2^) - {g{x) + 1} 



25{x) 



(45) 



where (?(x) = ln(2x) + 7. This expression for P2{x) is same as in Eq. ( l3Ti) . 
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The mass distribution obtained using perturbative approach upto 2"'^ order is 
compared to that of the Monte-Carlo simulations in Fig. |2l 

In the following section we study the model ASCM with random sequential 
dynamics. 

3.3. Random sequential dynamics : 

In this subsection we study the asymmetric sticky chipping model (ASCM) with 
random sequential update where at each time step a site i is chosen randomly, and from 
this site, the mass (1 — \)riXi is transported to the right neighbor i + 1. Rest of the 
mass {A + (1 — A)(l — rj)}xj remain at the departure site i. As usual, is a random 
number distributed uniformly in the interval (0, 1). 

Like the parallel dynamics, here too one can obtain all the moments exactly in 
terms of A assuming that the steady state has a product measure. Then we use the 
perturbative approach similar to the previous case of parallel dynamics. Clearly the 
moments (x") are given by 

(xD = i({A + (1 - A)(l - r.)rxn + \{{x, + (1 - A)r,_ix._i}"), (46) 

for alH = 1, 2, . . . , L, where, the 1** term represents that on the average half of the time 
Xxi amount of the mass sticks to the site i and the rest is chipped off. Of this chipped 
off mass (1 — A)xj, the amount (1 — X)riXi is transported to site i + 1 and the rest is 
retained at the departure site i. Thus half of the time site i retains {A + (l — A)(l — rj)}a;j 
amount of mass and in the other half i receives (1 — A)rj_iXi_i mass from site i — 1. 



3.3.1. Approach I : 



We begin with Eq. 

(a;") 
(n+1)! " 



n 



-1 



6]) and (x") in terms of lower moments {x'^) with m < n : 
n + l)A + A"+i-(l-A)"+^ 



n-l 



X Y.{n - k + - X) 

k=l 



fc+l 



(A; + 1)! {n-k + l)\' 



(47) 



After expanding the moments (x") as power series in A about A = as (x") 



oo \ fc/^ 
^fc 



in) 



where the coefficients C^"^ are independent of A. For convenience we 



choose C, 



in) 



(n + 1)\A 



in) 



f, — V"- -r ^j-^k ■ Using these in Eq. fll7|) and then equating the coeffi- 
cients of different powers of A, we obtain a set of recursive relations which are needed 



to be satisfied by A/^ 



for different orders. However obtaining the solution for 



those recursive relations, is usually difficult. Instead one can use Vk{s), the generating 
)r . 

^-^ms)] = ^-'[£{sv,{s)}]. 



functions for A^^^ (see Eq. fHT]) ) and then calculate Pfc(x) by using Eqs. f lT2|) and pUj) 

Pkix) 
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Qth Qj^dg^ perturbation in X : 



Coefficient of terms independent of A in Eq. ( H7|) gives 



n-l 



(k) An-k) 







(4J 



k=l 



Multiplying both sides of Eq. ( HHj) by (— s)" and taking sum over n we obtain an 
equation in terms of the generating function Vo(s) define in (l4Tll . 

sVo{s)V^{s) - 2syo'(s) + \/o2(s) _ v,{s) = 0. 

Since Vo(s) satiesfies the boundary condition (|39|1 . we get the solution, 



which gives 



V2^+l - 1 



d 



-x/2 



(49) 



(50) 



ds \r^^ 

Note that Pq{x) is same as the mass distribution of the ASCM with A = 0, obtained 
earlier [17] for random sequential dynamics. The steady state distribution function for 
ASCM with random sequential dynamics obtained from Monte-Carlo simulation of sys- 
tem size L = 100 for A = 0, 0.1 and 0.2, is shown in Fig. |3]^a). In Fig. Et^b) we compare 
Po{x) with P{x) obtained from Monte-Carlo simulations. 




Figure 3. (a) The steady state distribution P{x) of ASCM with random sequential 
update obtained from Monte-Carlo simulations of system size L = 100 for A = 0,0.1 
and 0.2. (b) P{x) for A = is compared with Eq. ([SC 



1^* order perturbation in A 



Collecting the coefficient of A in Eq. ( 1471) one can get 



(n - 1)A[ 



k=l 



{n-k + l){k + 1)4''^4"~''^ + in-k + l) 

k=l " 

.{n-k) 
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This recursion relation can be converted to a differential equation in terms of the 
generating functions Vo(s) and Vi{s). 

- sVo{s)V;is) + 2sV;is) - sV^{s)Vi{s) - 2Vo{s)Vi{s) + Vi{s) + sVo{s)V^{s) 
+ sVo'2(s)-K)(s) + l = 

where Vo(s) is given in Eq. (149|1 . Using the boundary condition (139|1 the solution is 
obtained as 

2s (s - 3) + 4 (s - 1) + 2 (s - s + 1) ln(s) 



where, s = a/1 + 2s. Thus we get 

d 



2ss 



(51) 



Pi{x) = 



-1 



ds 



72 



^ {-x^l\g{x) + 2} + x-^l\g{x) + 4}| + 4 



- 2b{x) 
(52) 

where g{x) = ln(2a;) + 7. In the above expression the term —26{x) is needed to ensure 
the condition dx Pi{x) = 0. 



The perturbative results up to order correction are then compared with P{x) 
obtained from the Monte-Carlo simulations of the model with different A. 




Figure 4. P{x) obtained from Monte-Carlo simulations of ASCM with random 
sequential dynamics on a lattice of size L — 100 are compared with the perturbative 
results for A = 0.1 and 0.2. In the upper ((a) and (b)) and lower ((c)and (d)) panel 
1** and 2"*^ order correction are being taken into account respectively. 



Fig. m shows that for A = 0.1, the perturbatively obtained result is very close to that 
obtained from Monte-Carlo simulations whereas for A = 0.2 there is a fair amount of 
discrepancy which may be improved by adding higher order corrections. 
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2nd Qj^d^j- perturbation in A ; 

From Eq. (H7|) collecting the coefficient of we have the equation 

n-l 



k=l 



,(fc) ^(n-fc) 



(k) ^(n-fe) 



(53) 



Using the generating functions Vo(s), V^i(s) and V2(s) from Eq. fHT]) . we get 



-is3yo'(^)V^,"(^) - |^5Vo(s)Fo"(s) + 2sVo'(s)V;'(s) + sVi{s)V^{s) - S\/2(^5)K)'(^) + 
syo(s)K'(5)-5K,(s)y2'(s)-2V^,(s)^(s)-Vb(s)-sV'i(s)K(s)-ri(s)2-Vl(s) + 2sy2'(s) + 

\/2(s)- 1 + 1 = 0. 

Since Vo(s) and Vi(s) are already known from Eq. fll^ and fl5T]) respectively, the solu- 
tion of the above equation can be obtained using the boundary conditions ( !39|l : 



Vo(s) 



l-S 



24s4(l+s) 



16 (s - 1) - 4s{s (4s + 9s + 6) + 26s - 6} 



+ 4{3 - 3s (s - 2) } In^(s) + 24{s (4s - 2) + s - 1} In s 
which results in 



P2(X) = ^ 



-1 



as 



-e 2 



8V27T 



--X 



-3/2 



l=(x^/^{Mx) + 2}^ + 6-^} 



16V27r 

^2 



-4xi/2{{(?(x) + 3} g 2 
-^x{(7(x) - 2 ln(2) + 3} + ^{(?(x) - 2 ln(2) + s} 



A5{x 



3 3 ^ ^ 



(54) 



where g{x) = ln(2x)+7. In the above expression the last two terms ensure normalization 
of probability and conservation of mass respectively (see Eq. (jl])). 

In Fig. m (c) and (d) we compare ^0(2^) + APi(x) + X'^P2{x) with the mass distri- 
bution obtained from Monte-Carlo simulations of a system of size L = 100. 



3.3.2. Approach II (limp^o) ■' 

As discussed in section 13.21 ASCM with random sequential dynamics can 
be considered as a special case (p — > limit) of parallel update scheme where the 
distribution of the random number is given by Eq. ( 132|) . In this limit terms iff{p'^) can 
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be neglected in Eq. (15^ . which gives 
d 
ds 



^{sVis)}[2-V{sil-\) 



{\s{l-r)Y 



'0 ^0 



ds 



(55) 



Solving this equation perturbatively we obtain the same Po{x), Pi{x) and ^2(2;) which 
turn out to be same as Eqs. fl50|) . fl52l) and fl5^ respectively. 

3.3.3. An analytically solvable case, A = 1/2 : 

It is interesting to note that the asymmetric sticky chipping model (ASCM) with 
continuous time update, can be solved analytically at A = 1/2. In this subsection we 
study this case in details. The dynamics is as follows. At each time step one site i with 
mass Xi is chosen randomly. Half of the amount (y) sticks to the site i. Of the rest 
half, is transported to the right neighbor i + 1 and rest is retained at site i. In 
other words the mass yrj is transported to site (i + 1) from site i and the mass 
remain at the departure site i. In this case we obtain P{x) as follows. From Eq. (H^ 
we get the equation 



2\[ 2 J / 2 \l 2 

Simplifying we get an equation that expresses the n^^ term in terms of lower moments 
up to n — 1: 

(^ + 1)! 2^= (fc + l)!(n-fc + l)!' ^ ' 

With the boundary condition (x) = 1 the solution is obtained as 

2n 



- ^ (57) 



and hence 

P(x) = ^-^y =4xe-2-. (58) 

^0 ^■ 

Thus in this case we obtain the closed form expression for the steady state mass 
distribution function. To check whether the result is consistent we compare our result 
with Monte-Carlo simulations of a system of size L = 100. 

Fig. dM*^)) shows that Eq. fISS]) agrees quite well with the distribution P{x) 
obtained from Monte-Carlo simulations. We have checked that the small discrepancy 
near the peak (shown in the inset of Fig. ([5])) does not disappear when the system size 
is increased. This small discrepancy observed indicates that possibly the steady state 
in this case is not actually factorized. 
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Figure 5. (a) The steady state distribution function P{x) (dashed line) of ASCM 
with random sequential dynamics for A = 1/2 obtained form Monte-Carlo simulations 
of system size L = 100, is compared with Eq. ([55]) (solid line). 



4. Symmetric Sticky Chipping Model (SSCM) : 

In this section we study the symmetric version of the model, namely SSCM. Here, 
the dynamics is as follows. In each time step at site i with mass Xi, A fraction of the 
mass sticks to the site and (1 — A)xj is chipped off and distributed randomly among 
both the neighbors of i i.e. {i ± 1) in a way that one of those receives rj(l — X)xi and 
the other receives the rest, (1 — Tj)(l — A)xj, where rj is a random number uniformly 
distributed in the interval (0,1). We will use two different update schemes - parallel 
update and random sequential update. 

4.1. Parallel Update : 

In this subsection we study the model where all the sites are updated parallelly 
(synchronously) using the dynamics mentioned above. Mathematically, the dynamics 
can be explicitly written as 



where, the 1** term on the r.h.s. represents the mass that sticks to the site i during the 
update, 2"''^ and 3'"'^ term correspond to the mass which site i receives from its neighbors 
i =F 1 respectively. 

The steady state mass distribution function can be obtained in the same approach as 
shown in case of ASCM. Here we follow the approach as shown in section [3^ Converting 
the dynamical equation for the system (Eq. fl59|) ) into an integral equation of P{x) (as 
shown in case of ASCM) and taking the Laplace transform of both sides of that equation, 
we get 



Xi{t + 1) = \xi{t) + (1 - A)[rj_iXi_i(t) + ri+iXi+i{t)] 



(59) 



Q{s) = QiXs)v\il - X)s) 




(60) 



- s=Xs 
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Note that V{s) is also a function of A. Here we have dropped the argument A in V{s) for 
notational convenience. We obtain the solution of this equation using the perturbative 
approach by taking the Taylor series expansion of both sides of Eq. ( l60l) . and equating 
the coefficients of different powers of A order by order. 



order perturbation in A : 

From Eq. ( l60l) we equate the coefficient of the terms independent of A, which 
gives 

- Voisf + Vo{s) = 
This can be solved using the boundary conditions ( l39i) : 

Vo{s) ^ 

This results in 



s + 2 



Axe 



-2x 



This expression is identical to that obtained in case of asymmetric sticky chipping model 
with parallel update at A = and the expression obtained in case of asymmetric sticky 
chipping model with random sequential update at A = 1/2. The steady state distribution 
finction P{x) obtained from Monte-Carlo simulations of SSCM with parallel dynamics 
for A = 0,0.1 and 0.2 with system size L = 100 is shown in Fig. [6](a). Fig. [6](b) 
corresponds to the comparison of Po{x) with Monte-Carlo simulations. 
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Figure 6. (a) shows P{x) obtained from Monte-Carlo simulations of SSCM with 
parallel dynamics for system size L — 100 for A = 0,0.1 and 0.2. (b) compares Po{x) 
with P{x) obtained from Monte-Carlo simulations. 



1 order perturbation in A .■ 



In Eq. ( l60l) collecting the coefficients of A^ we get the differential equation 



s{2 + sfv;{s) + (2 + s)\s - 2)Vi{s) + 4s^ 
which can be solved using the boundary conditions fl39l) : 

4s In(^) 







[s + 2) 
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Thus we get 



-1 



16e 



-2x 



where g{x) = ln(2x) + 7. We compare the perturbatively obtained distributions (sohd 
hues) with 1** and 2"*^ order correction, with P{x) obtained from Monte-Carlo simula- 
tions of system size L = 100, for A = 0.1 in Fig. [Tl^a) and (b) respectively. 





4 5 



Figure 7. The distribution functions obtained perturbatively (solid lines) with 1** 
and 2""^ order correction are compared with the results obtained from Monte-Carlo 
simulations (dashed lines) of SSCM with parallel dynamics of system size L = 100 for 
A = 0.1. 



2nd Qj^fi^j. perturbation in A .■ 



A differential equation in V2{s) can be obtained from coefficient of A in Eq. ([60 

s{2 + sYVi{s) + (2 + sf{s - 2)V2{s) - 16s In (^) 

+ s In (i±^) + 2} - s2{s(3s + 4) + 24} = 0. 
Using the boundary conditions fl39|) one can obtain the solution. 



9s^ + (36 - 47r2)s - Svr^ + 241n(f){51n(-f) + 2 - 5} + 12(s - 2) In^d 



+ 24SLi2{{ 

where, s = s + 2, Li2{z) = 
Thus we obtain 

— Q„-2x 



E 



m=l rr)2 



, is the Polylog function [201 HI 



Ax^{g\x) - 2g{x) + Ah{2x) + + Sx^j - g\x) + \g{x){2 - Ei{2x)} 

- 2h{2x) + i/(2a;) - | + } + 2x{g\x) + 2g{x){l + Ei{2x)} - ^ - 2f{2x) 
+ e2^Ei(-2x)} - {g{x) + 1} - e^''Ei{-2x) + 65(x) 



where, g{x) = ln(2a;) + 7, Ei{x) = —J^^dt is the exponential integral [20] . 
h{x) = E^=o^^;t^ and fix) = E^^.^b + ^'^'Km)]^ with ^(o)(m) = ^j^. 
In the expression for P2{x) the term 6S{x) arises in order to satisfy the normal- 
ization condition ([!]). The mass distribution function up to 2""^ order correction 
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P{x) = Po{x) + XPi{x) + X'^P2{x) is compared with P{x) obtained from Monte-Carlo 
simulations of a system of size L = 100, in Fig. Wljo)- The figure shows that compared to 
the 1** order correction, the perturbative result up to the 2"*^ order correction is closer 
to P{x) obtained from Monte-Carlo simulations. For better accuracy one has to go to 
higher orders. 

In the following section we study SSCM with random sequential dynamics. 



4-2. Random sequential update : 

In this subsection we study the symmetric sticky chipping model using random 
sequential update. The dynamics is as follows. In each step choose one site i randomly. 
Of the mass Xj, the amount Axj sticks to the site, (1 — A)rjXj is transported to the right 
neighbor i + 1 and (1 — A)(l — rj)xj is transported to the left neighbor i — 1. Assuming 
that the steady state has a product measure, the moments (x") can be expressed as 

(x") = i((Ax,)") + i(l - A)" ((x, + r,_ix,_i)" + {x, + (1 - n+i)x,+i}") 

= ^((Ax.)") + ^(1 - A)"((x. + r,_ix._i)") (61) 

i.e., on the average one third of the time Axj sticks to the site i, one third of the 
time (1 — A)rj_iXj_i is transported to site i from site i — 1 and one third of the time 
(1 — A)(l — ri+i)xi+i is transported to site i from site i + 1. 
Now, consider the case A = 0. Then, Eq. (!6T]) takes the form 

(x") = ^((x, + r,_ix,_i)"). 

Starting from this equation the higher moments can be expressed in terms of the lower 
moments in the form of recursive relation. 

2 (n + 1)! ^^^"^ ' {n-k + l)\{k + l)\ 

or, 

n-l 

{n - 1)4") = 2 ^ (n - + 1)4')4"-') (62) 

k=l 

where, A'^^ = j^^- The generating function for A*^^ is Vo(s) = Z]J^Lo("~'^)"^o"^- 
Multiplying both sides of Eq. ( 1621) by {—s)^ and taking sum over n we obtain a relation 
in terms of the generating function Vo(s): 
_ l-3Vois) + 2V,'is) 
s{3-2Vo{s)} ■ 

This can be solved using the boundary condition Eq. ( 15^ and we get 

Ms) = 

which results in 

' d 



Poix) = 



^^{sVois)} 



x-y^e-'^/^ + Six). (63) 
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In the above expression the term 6{x) is needed to ensure the normahzation condition 

We compare the perturbative result with the Monte-Carlo simulations of a system 
of size L = 100 in Fig. E 



0.8 

^ 0.6 
0.4 
0.2 


Figure 8. Comparison of Eq. (solid line) with P{x) obtained from Monte-Carlo 
simulations (dashed line) of SSCM with random sequential dynamics for A = of a 
system of size L = 100. 

The Fig. |H] shows that Eq. (16 3 p does not agree with the result obtained from 
Monte-Carlo simulations. Thus, in this case the assumption that the steady state has a 
product measure, fails badly. We do not proceed further in this case. 




5. Summary and conclusion : 

In this article we have studied conserved mass transfer processes using different 
models in one dimension with the mass treated as a continuous variable. In these 
chipping models at each step of update the mass at a site is distributed randomly 
among the site and its neighbors. We introduced a stickiness parameter A (0 < A < 1) 
denoting the fraction of mass that sticks to the site. Thus, a fixed fraction (1 — A) 
is transported to the neighbors. We have studied the models with asymmetric and 
symmetric dynamics with both parallel and random sequential update. 

Finding the stationary state distribution analytically, is not always possible. We 
introduce a perturbative approach to obtain the approximate distribution function 
assuming that the steady state has a product measure. We expand the distribution 
function about A = as P{x) = ^2^=0 ^''^kix), and obtain the distributions up to 2"'^ 
order correction, P{x) = Po{x) + APi(x) + \'^P2{x) . This gives the correct values of all 
the moments ((x")) up to the order A^. We compare our perturbative results with that 
obtained from Monte-Carlo simulations, which shows that the perturbative results agrees 
quite well with that obtained from Monte-Carlo simulation. This indicate that product 
measure assumptions for the steady state holds well. For the special case of A = in 
all the models studied in this article and for the asymmetric sticky chipping model with 
random sequential update at A = 1/2 , we have obtained the closed form expression of 
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the distribution function. However, Fig. clearly indicate that the assumption that 
the steady state has a product measure, fails badly in case of symmetric sticky chipping 
model with random sequential update. 

From the perturbative series of P{x) in A up to 2"^^ order it is not always possible 
to find a general expression for Pk{x). In the sticky chipping models studied in 
this article, the mass distribution can be expressed generically in the form Pk{x) = 
Sz,m G'fc'"''^' In™" (x) where G^.'™''' are constant coefficients. Usually the sum runs over 
finite number of terms (here I < k + 1 and m < k). Thus for small k, one can always 
obtain G]^"^^ explicitly. Getting a generic expression, however, is not always possible. 
The existences of the logarithmic term in the series expansion, illustrates the complexity 
and "slow" convergence of P{x) to a limiting function. 
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